Although in the modern world dieting is a hot and widespread topic due to the implications it has on health (obesity pandemic) and aesthetics (especially as social media constantly convey unrealistic expectations around body image), this blog post will reduce its discussion and analysis to the mathematical formulations that can arise from the challenge of choosing what and how much to eat given different criteria/objectives.
My first encounter with this problem was in Boyd's book (pg. 162), where the objective is to build a diet that contains enough nutrients to be a healthy one but is the cheapest possible. What does the "but" in the last sentence for? Because the cheapest diet without restrictions is to eat nothing 🤡.
For this study, the dataset used comes from this page, and the daily intake can be found in ./data/daily_intake.csv on the code repository.
You can take a look at the data on the following cell:
import pandas as pd
df = pd.read_csv('./data/stigler.csv')
print("df:")
display(df)
daily_intake = pd.read_csv('./data/daily_intake.csv')
print("daily_intake:")
display(daily_intake)
df:
| commodity | unit | price_cents | calories_k | protein_g | calcium_g | iron_mg | vitamin_a_kiu | vitamin_b1_mg | vitamin_b2_mg | niacin_mg | vitamin_c_mg | |
|---|---|---|---|---|---|---|---|---|---|---|---|---|
| 0 | Wheat Flour (Enriched) | 10 lb. | 36.0 | 44.7 | 1411 | 2.0 | 365 | 0.0 | 55.4 | 33.3 | 441 | 0 |
| 1 | Macaroni | 1 lb. | 14.1 | 11.6 | 418 | 0.7 | 54 | 0.0 | 3.2 | 1.9 | 68 | 0 |
| 2 | Wheat Cereal (Enriched) | 28 oz. | 24.2 | 11.8 | 377 | 14.4 | 175 | 0.0 | 14.4 | 8.8 | 114 | 0 |
| 3 | Corn Flakes | 8 oz. | 7.1 | 11.4 | 252 | 0.1 | 56 | 0.0 | 13.5 | 2.3 | 68 | 0 |
| 4 | Corn Meal | 1 lb. | 4.6 | 36.0 | 897 | 1.7 | 99 | 30.9 | 17.4 | 7.9 | 106 | 0 |
| ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... |
| 72 | Chocolate | 8 oz. | 16.2 | 8.0 | 77 | 1.3 | 39 | 0.0 | 0.9 | 3.4 | 14 | 0 |
| 73 | Sugar | 10 lb. | 51.7 | 34.9 | 0 | 0.0 | 0 | 0.0 | 0.0 | 0.0 | 0 | 0 |
| 74 | Corn Syrup | 24 oz. | 13.7 | 14.7 | 0 | 0.5 | 74 | 0.0 | 0.0 | 0.0 | 5 | 0 |
| 75 | Molasses | 18 oz. | 13.6 | 9.0 | 0 | 10.3 | 244 | 0.0 | 1.9 | 7.5 | 146 | 0 |
| 76 | Strawberry Preserves | 1 lb. | 20.5 | 6.4 | 11 | 0.4 | 7 | 0.2 | 0.2 | 0.4 | 3 | 0 |
77 rows × 12 columns
daily_intake:
| nutrient | daily_recommended_intake | unit | |
|---|---|---|---|
| 0 | Calories | 3.0 | k Calories |
| 1 | Protein | 70.0 | g |
| 2 | Calcium | 0.8 | g |
| 3 | Iron | 12.0 | mg |
| 4 | Vitamin A | 5.0 | KIU |
| 5 | Vitamin B1 | 1.8 | mg |
| 6 | Vitamin B2 | 2.7 | mg |
| 7 | Niacin | 18.0 | mg |
| 8 | Vitamin C | 75.0 | mg |
As data preprocessing steps we took df and:
If you want to see the code go to the GitHub repository.
from dietopt.preprocess import apply_preprocess
df_processed = apply_preprocess(df)
print("df_processed:")
df_processed
df_processed:
| commodity | unit | USD | calories_k | protein_g | calcium_g | iron_mg | vitamin_a_kiu | vitamin_b1_mg | vitamin_b2_mg | niacin_mg | vitamin_c_mg | |
|---|---|---|---|---|---|---|---|---|---|---|---|---|
| 0 | Wheat Flour (Enriched) | lb. | 7.718292 | 4.470000 | 141.100000 | 0.200000 | 36.500000 | 0.0 | 5.540000 | 3.330000 | 44.100000 | 0.0 |
| 1 | Macaroni | lb. | 3.022998 | 11.600000 | 418.000000 | 0.700000 | 54.000000 | 0.0 | 3.200000 | 1.900000 | 68.000000 | 0.0 |
| 2 | Wheat Cereal (Enriched) | oz. | 5.188407 | 0.421429 | 13.464286 | 0.514286 | 6.250000 | 0.0 | 0.514286 | 0.314286 | 4.071429 | 0.0 |
| 3 | Corn Flakes | oz. | 1.522219 | 1.425000 | 31.500000 | 0.012500 | 7.000000 | 0.0 | 1.687500 | 0.287500 | 8.500000 | 0.0 |
| 4 | Corn Meal | lb. | 0.986226 | 36.000000 | 897.000000 | 1.700000 | 99.000000 | 30.9 | 17.400000 | 7.900000 | 106.000000 | 0.0 |
| ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... |
| 72 | Chocolate | oz. | 3.473231 | 1.000000 | 9.625000 | 0.162500 | 4.875000 | 0.0 | 0.112500 | 0.425000 | 1.750000 | 0.0 |
| 73 | Sugar | lb. | 11.084325 | 3.490000 | 0.000000 | 0.000000 | 0.000000 | 0.0 | 0.000000 | 0.000000 | 0.000000 | 0.0 |
| 74 | Corn Syrup | oz. | 2.937239 | 0.612500 | 0.000000 | 0.020833 | 3.083333 | 0.0 | 0.000000 | 0.000000 | 0.208333 | 0.0 |
| 75 | Molasses | oz. | 2.915799 | 0.500000 | 0.000000 | 0.572222 | 13.555556 | 0.0 | 0.105556 | 0.416667 | 8.111111 | 0.0 |
| 76 | Strawberry Preserves | lb. | 4.395138 | 6.400000 | 11.000000 | 0.400000 | 7.000000 | 0.2 | 0.200000 | 0.400000 | 3.000000 | 0.0 |
77 rows × 12 columns
As now we have the data, the mathematical formulation may be clearer.
Our diet must contain $m$ nutrients with at least $b_1$, . . . , $b_m$ quantities, as you may see
on the daily_intake dataframe we have $m=9$ and our vector $b$ will be the values in column
daily_recommended_intake of daily_intake dataframe. To build the diet we choose nonnegative quantities $x_1$, . . . , $x_n$ of
$n$ different foods, in our case $n=77$ which is the number of rows of our df_processed dataset.
One unit quantity of food $j$ contains an amount $a_{ij}$ of nutrient
$i$ and has a cost of $c_j$ (note that the values in df_processed are the transpose of the following matrix $A$).
We want to determine the cheapest diet, $x_{opt}$, that satisfies the nutritional requirements, daily_recommended_intake. This problem
can be formulated as the linear program (LP):
Where if $a, b \in \R^n_+$, $a \succeq b$ means that $\forall i = 1, . . . , n \implies a_i \ge b_i$.
Now that we have the data and the mathematical formulation, we can model the problem and find the "optimal" diet given our criteria (the objective function). To do so, we will use CVXPY, a python package for solving convex problems. Although our problem is a linear one and CVXPY maybe feel like overkill, as linear problems are convex ones, this library is a more general tool than using a linear solver, meaning that we can keep using the same tool for other non-linear problems in the future.
The following code uses CVXPY to solve the problem of Eq. 1:
import cvxpy as cp
# Problem data
c = df_processed.USD.values
A = df_processed.iloc[:, 3:].T.values
b = daily_intake.daily_recommended_intake.values
# Constructing the problem
x = cp.Variable(shape=len(c))
objective = cp.Minimize(c@x)
constraints = [
A@x >= b,
x >= 0,
]
problem = cp.Problem(objective=objective, constraints=constraints)
result = problem.solve()
print(f"{result = }")
print(f"{x.value = }")
result = 0.1536070376685601
x.value = array([-3.76573316e-14, 8.68348920e-14, -1.46738328e-13, 8.13624289e-14,
5.34424634e-03, 2.95446061e-14, 3.06770764e-13, 4.37997030e-13,
2.77853882e-13, 2.04462612e-13, 1.77817196e-13, -1.69845507e-14,
6.44932975e-14, 3.49824397e-13, 9.76998354e-14, -7.88504668e-14,
2.33506059e-14, -1.23905986e-14, 4.03105827e-14, 4.02896389e-14,
1.10932580e-13, -1.13920983e-14, -4.78824654e-14, 2.04427572e-13,
-3.59915089e-14, -2.24149401e-14, -1.19513449e-14, 2.23277626e-14,
7.37335280e-14, 5.53354150e-14, -1.12987895e-15, -3.96683084e-14,
-6.50090199e-15, 1.99263707e-14, -1.15558487e-14, -2.94876516e-15,
4.50554631e-14, -1.87698517e-14, -4.19150510e-14, -2.37919029e-14,
4.57123757e-13, 2.85788206e-13, 2.81256178e-14, 3.21319703e-14,
2.43133822e-13, 1.13132451e-02, 4.68001567e-13, 1.77437607e-13,
1.82422902e-13, 8.41262040e-13, -4.33849115e-14, 5.17574850e-03,
2.39843674e-12, 2.17471565e-14, -4.93345153e-14, 2.21769476e-14,
-1.03700236e-14, 1.33085121e-13, 6.30155777e-14, 1.06124977e-13,
7.10903546e-14, 1.88576706e-13, 1.03854535e-13, 5.39432412e-14,
1.95543157e-13, -5.32264168e-14, 3.87732822e-13, 2.89602726e-13,
1.03066891e-01, -3.54310802e-14, 4.52763432e-14, 3.49870292e-14,
-1.32270420e-13, -1.99661357e-13, -1.94840287e-13, -8.99133127e-14,
-1.36720909e-13])
Translating the numbers, or the $x$ vector, of our problem, the optimal diet (weekly to make the numbers more readable) would be:
from dietopt.utils import compute_weekly_diet_df
df_weekly_diet = compute_weekly_diet_df(df=df_processed, qtds=x.value.tolist())
df_weekly_diet
| commodity | qtd_weekly | usd_weekly | |
|---|---|---|---|
| 0 | Navy Beans, Dried | 0.7215 lb. | 0.91 |
| 1 | Cabbage | 0.0792 lb. | 0.06 |
| 2 | Spinach | 0.0362 lb. | 0.06 |
| 3 | Corn Meal | 0.0374 lb. | 0.04 |
| 4 | TOTAL | 1.07 |
Ok ok, a diet that gives us enough nutrients and costs 1.07 USD/week seems pretty cheap, as was stated from the optimization problem. But, looking at the menu... will I just eat 4 different items? Yeah... this does not seem like a diverse diet that I'll be able to stick to. Can I model something like diversity on my diet?
x = cp.Variable(shape=len(c))
objective = cp.Minimize(c@x + 10*cp.norm(x,2)**2)
constraints = [x >= 0, A@x >= b]
problem = cp.Problem(objective=objective, constraints=constraints)
result = problem.solve()
df_weekly_diet = compute_weekly_diet_df(df=df_processed, qtds=x.value.tolist())
df_weekly_diet
| commodity | qtd_weekly | usd_weekly | |
|---|---|---|---|
| 0 | Navy Beans, Dried | 0.3731 lb. | 0.47 |
| 1 | Lima Beans, Dried | 0.1857 lb. | 0.35 |
| 2 | Corn Meal | 0.1935 lb. | 0.19 |
| 3 | Cabbage | 0.0785 lb. | 0.06 |
| 4 | Spinach | 0.0234 lb. | 0.04 |
| 5 | Peas, Dried | 0.0169 lb. | 0.03 |
| 6 | Sweet Potatoes | 0.0205 lb. | 0.02 |
| 7 | TOTAL | 1.16 |
The little trick applied above was to add norm 2 regularization to the $x$ variable, the new optimization problem would be:
$$ \begin{equation} \tag{2} \begin{split} \text{minimize } \ & c^tx + 10\lVert x \rVert _2^2\\ \text{subject to } \ & Ax \succeq b \\ & x \succeq 0 \end{split} \end{equation} $$But why the $l_2$-norm was applied? In a nutshell, if you derive the $\lVert x \rVert _2^2$ by some component of $x$, like $x_i$, you will find
$$ \frac{\partial}{\partial x_i} \lVert x \rVert _2^2 = 2x_i $$Meaning that the bigger is $x_i$ biggest is the incentive (derivative) to decrease its value on the objective function. As $x_i$ gets close to zero the incentive decreases too, meaning, in practical terms, that we penalize big (relative) values and don't care too much about small ones. That would be something like: "I want a cheap diet but don't want to just eat one thing, as I want some variety, even if I can eat just a little".
Note: where does this random $10$ come from? It was chosen by try-and-error until a found a more diverse diet that didn't cost too much (find the balance between both objectives).
An interesting thing happened by the end of the previous section: we noted that our objective wasn't really aligned with what we wanted. This is an important conclusion that you might not see much value in, but we real-life we often see this misalignment between the mathematical formulation and the business need. When this misalignment happens, solving the optimization problem does not translate into solving the real-world problem.
So, initially, we wanted a diet that had all the necessary nutrients while being cheap as possible (see Eq. 1), but, what if, we were looking for getting lean, what would happen?
This time we want to minimize our calories while ingesting enough of the other nutrients. The $\sigma$ vector represents how many calories each food contains.
$$ \begin{equation} \tag{3} \begin{split} \text{minimize } \ & \sigma^tx \\ \text{subject to } \ & Ax \succeq b_{\sigma} \\ & x \succeq 0 \end{split} \end{equation} $$On equation 3, $b_{\sigma}$ is the modified version of our daily_intake where the minimum value for calories is $0$ (zero).
import numpy as np
A = df_processed.iloc[:, 3:].T.values
b = daily_intake.daily_recommended_intake.values.copy()
b[0] = 0 # no need for calories constraint
calories = df_processed.calories_k.values
x = cp.Variable(shape=len(calories))
objective = cp.Minimize(calories@x)
constraints = [x >= 0, A@x >= b]
problem = cp.Problem(objective=objective, constraints=constraints)
result = problem.solve()
print(f"Calories: {result:.2f} kcal/day")
cost = df_processed.USD.values
df_weekly_diet = compute_weekly_diet_df(df=df_processed, qtds=x.value.tolist())
df_weekly_diet
Calories: 0.60 kcal/day
| commodity | qtd_weekly | usd_weekly | |
|---|---|---|---|
| 0 | Salmon, Pink (can) | 6.5861 oz. | 18.36 |
| 1 | Coffee | 3.3034 lb. | 15.86 |
| 2 | Liver (Beef) | 0.4617 lb. | 2.65 |
| 3 | Tea | 0.4397 lb. | 1.64 |
| 4 | Celery | 0.9029 stalk | 1.41 |
| 5 | TOTAL | 39.92 |
It would be a more expensive diet but with a really small quantity of calories daily, only 600 cal.
Disclaimer: Please, remember this is a mathematical exercise, not a real diet, don't be dumb, and just eat salmon from now on.
Thinking about dieting we may realize that one of the big challenges of it is that people need to change their eating habits, meaning, they end up eating what they don't like. So, let's add a new column to our data set with the "tastiness" of the commodities, ranging from -5, "I really despise it", to +5, "If I could, I would eat this every day".
To avoid eating too much we also added an upper bound to the nutrients in the diet (which include calories, probably a big effect by my tastiness preferences 😁). The $\tau$ letter will represent the vector of the tastiness of each commodity.
$$ \begin{equation} \tag{4} \begin{split} \text{minimize } \ & -\tau^tx \\ \text{subject to } \ & Ax \succeq b \\ & Ax \preceq 1.2b \\ & x \succeq 0 \end{split} \end{equation} $$PS: The negative signal on the objective function of Eq. 4 is to translate the tastiness maximization problem to a minimization one.
On the next cell, you may see the tastiness I see on the commodities available:
tastiness_df = pd.read_csv('./data/tastiness_by_romulo.csv')
tastiness_df
| tastiness | commoditie | |
|---|---|---|
| 0 | 0 | Wheat Flour (Enriched) |
| 1 | 3 | Macaroni |
| 2 | 2 | Wheat Cereal (Enriched) |
| 3 | 3 | Corn Flakes |
| 4 | 0 | Corn Meal |
| ... | ... | ... |
| 72 | 5 | Chocolate |
| 73 | 0 | Sugar |
| 74 | -2 | Corn Syrup |
| 75 | 3 | Molasses |
| 76 | 3 | "Strawberry Preserves" |
77 rows × 2 columns
Now, the code for solving the problem:
A = df_processed.iloc[:, 3:].T.values
b = daily_intake.daily_recommended_intake.values.copy()
tastiness = tastiness_df['tastiness'].values
x = cp.Variable(shape=len(calories))
objective = cp.Minimize(-tastiness@x)
constraints = [
x >= 0,
A@x >= b,
A@x <= 1.2*b
]
problem = cp.Problem(objective=objective, constraints=constraints)
result = problem.solve()
print(f"{result = }")
calories = df_processed.calories_k.values
print(f"Calories: {calories@x.value:.2f} k/day")
df_weekly_diet = compute_weekly_diet_df(df=df_processed, qtds=x.value.tolist())
df_weekly_diet
result = -12.753855640167883 Calories: 3.60 k/day Cost: 71.15 USD/week
| commodity | qtd_weekly | usd_weekly | |
|---|---|---|---|
| 0 | Chocolate | 13.4992 oz. | 46.89 |
| 1 | Pineapple (can) | 1.2057 No. 2 1/2 | 5.51 |
| 2 | Evaporated Milk (can) | 3.1825 oz. | 4.57 |
| 3 | Leg of Lamb | 0.6486 lb. | 3.84 |
| 4 | Salmon, Pink (can) | 1.3280 oz. | 3.7 |
| 5 | Coffee | 0.6892 lb. | 3.31 |
| 6 | Butter | 0.3695 lb. | 2.44 |
| 7 | Pork Chops | 0.1338 lb. | 0.88 |
| 8 | Spinach | 0.0112 lb. | 0.02 |
| 9 | TOTAL | 71.16 |
Wooww, 71.16 USD/week is kind of expensive, couldn't we add the cost to the objective function so we could control it a bit?
The new problem formulation would be:
$$ \begin{equation} \tag{5} \begin{split} \text{minimize } \ & -\tau^tx + c^tx\\ \text{subject to } \ & Ax \succeq b \\ & Ax \preceq 1.2b \\ & x \succeq 0 \end{split} \end{equation} $$A = df_processed.iloc[:, 3:].T.values
b = daily_intake.daily_recommended_intake.values.copy()
c = df_processed.USD.values
x = cp.Variable(shape=len(calories))
objective = cp.Minimize(-tastiness@x + c@x )
constraints = [x >= 0, A@x >= b, A@x <= 1.2*b]
problem = cp.Problem(objective=objective, constraints=constraints)
result = problem.solve()
print(f"{result = }")
calories = df_processed.calories_k.values
print(f"Calories: {calories@x.value:.2f} k/day")
df_weekly_diet = compute_weekly_diet_df(df=df_processed, qtds=x.value.tolist())
df_weekly_diet
result = -3.88875776195683 Calories: 3.07 k/day Cost: 44.21 USD/week
| commodity | qtd_weekly | usd_weekly | |
|---|---|---|---|
| 0 | Chocolate | 7.0156 oz. | 24.37 |
| 1 | Corn Flakes | 7.0783 oz. | 10.77 |
| 2 | Evaporated Milk (can) | 4.9459 oz. | 7.1 |
| 3 | Tea | 0.3031 lb. | 1.13 |
| 4 | Corn (can) | 0.2390 No. 2 | 0.53 |
| 5 | Liver (Beef) | 0.0367 lb. | 0.21 |
| 6 | Cabbage | 0.0714 lb. | 0.06 |
| 7 | Spinach | 0.0180 lb. | 0.03 |
| 8 | TOTAL | 44.2 |
Ok ok, 44.21 USD/week seems more reasonable. As an effect of adding the cost to our objective function we got rid of Pineapple although I like it (+5 in tastiness) it seems a bit expensive compared to others. Because we now consider cost, Liver also got into my diet although I'm not such a big fan (-2 in tastiness).
Now we have an interesting problem: as we have two objectives in the same objective function, they might compete. Do they compete? We can test it by optimizing for each objective separately:
A = df_processed.iloc[:, 3:].T.values
b = daily_intake.daily_recommended_intake.values.copy()
c = df_processed.USD.values
x = cp.Variable(shape=len(calories))
constraints = [x >= 0, A@x >= b, A@x <= 1.2*b]
objectives_dict = {
"Cost": cp.Minimize(c@x),
"Tastiness": cp.Minimize(-tastiness@x),
}
for tag, objective_function in objectives_dict.items():
print(f"Solving for: {tag}")
problem = cp.Problem(objective=objective_function, constraints=constraints)
result = problem.solve()
print(f"{result = }")
calories = df_processed.calories_k.values
print(f"Calories: {calories@x.value:.2f} k/day")
cost = df_processed.USD.values
print(f"Cost: {cost@x.value*7:.2f} USD/week")
df_weekly_diet = compute_weekly_diet_df(df=df_processed, qtds=x.value.tolist())
display(df_weekly_diet)
print('\n')
Solving for: Cost result = 0.3957707708368031 Calories: 3.00 k/day Cost: 2.77 USD/week
| commodity | qtd_weekly | usd_weekly | |
|---|---|---|---|
| 0 | Liver (Beef) | 0.1466 lb. | 0.84 |
| 1 | Milk | 0.3343 qt. | 0.79 |
| 2 | White Bread (Enriched) | 0.2973 lb. | 0.5 |
| 3 | Corn Meal | 0.2758 lb. | 0.27 |
| 4 | Lard | 0.0712 lb. | 0.15 |
| 5 | Onions | 0.1580 lb. | 0.12 |
| 6 | Cabbage | 0.0571 lb. | 0.05 |
| 7 | Peanut Butter | 0.0112 lb. | 0.04 |
| 8 | TOTAL | 2.76 |
Solving for: Tastiness result = -12.753855640167883 Calories: 3.60 k/day Cost: 71.15 USD/week
| commodity | qtd_weekly | usd_weekly | |
|---|---|---|---|
| 0 | Chocolate | 13.4992 oz. | 46.89 |
| 1 | Pineapple (can) | 1.2057 No. 2 1/2 | 5.51 |
| 2 | Evaporated Milk (can) | 3.1825 oz. | 4.57 |
| 3 | Leg of Lamb | 0.6486 lb. | 3.84 |
| 4 | Salmon, Pink (can) | 1.3280 oz. | 3.7 |
| 5 | Coffee | 0.6892 lb. | 3.31 |
| 6 | Butter | 0.3695 lb. | 2.44 |
| 7 | Pork Chops | 0.1338 lb. | 0.88 |
| 8 | Spinach | 0.0112 lb. | 0.02 |
| 9 | TOTAL | 71.16 |
As you can see, we have different results (different $x_{opt}$) depending on the objective. We can visualize it in a 2D plane where each axis is one objective:
import plotly.express as px
A = df_processed.iloc[:, 3:].T.values
b = daily_intake.daily_recommended_intake.values.copy()
c = df_processed.USD.values
x = cp.Variable(shape=len(calories))
constraints = [x >= 0, A@x >= b, A@x <= 1.2*b]
objectives_dict = {
"Cost": cp.Minimize(c@x),
"Tastiness": cp.Minimize(-tastiness@x),
}
x_list = []
y_list = []
for tag, objective_function in objectives_dict.items():
problem = cp.Problem(objective=objective_function, constraints=constraints)
result = problem.solve()
cost = c@x.value
tastiness_ = -tastiness@x.value
x_list.append(cost)
y_list.append(tastiness_)
fig = px.scatter(x=x_list, y=y_list)
size = 600
fig.update_layout(
xaxis_title='Cost',
yaxis_title='- Tastiness',
width=size, height=size
)
fig.update_traces(marker=dict(size=20))
fig.show()
As you may see, we don't have a "perfect solution" where there the left-most point is the downwards-most point too.
To analyze the trade-off between the objectives we may give different weights for each objective and explore the curve that lies between the graphs of the previous plot. Let's do this!
For finding the Pareto curve we will weigh the objectives using the variable $\mu$, as you may see in Eq. 6.
$$ \begin{equation} \tag{6} \begin{split} \text{minimize } \ & -(1-\mu)\tau^tx + \mu c^tx\\ \text{subject to } \ & Ax \succeq b \\ & Ax \preceq 1.2b \\ & x \succeq 0 \end{split} \end{equation} $$The variable $\mu \in \R$ will assume many values between $[0,1]$ so we can see the trade-off between the two objectives.
from tqdm import tqdm
A = df_processed.iloc[:, 3:].T.values
b = daily_intake.daily_recommended_intake.values.copy()
c = df_processed.USD.values
mu_list = [0] + np.logspace(-5, 0, num=1_000).tolist()
x = cp.Variable(shape=len(calories))
constraints = [x >= 0, A@x >= b, A@x <= 1.2*b]
cost_list = []
tastiness_list = []
for mu in tqdm(mu_list):
objective_function = cp.Minimize(-(1 - mu)*tastiness@x + mu*c@x )
problem = cp.Problem(objective=objective_function, constraints=constraints)
result = problem.solve()
cost_ = c@x.value
tastiness_ = -tastiness@x.value
cost_list.append(cost_)
tastiness_list.append(tastiness_)
100%|██████████| 1001/1001 [00:05<00:00, 197.82it/s]
import plotly.graph_objects as go
import plotly.io as pio
pio.renderers.default = 'plotly_mimetype+notebook'
fig = go.Figure(go.Scatter(x=cost_list, y=tastiness_list, mode='lines'))
size = 600
fig.update_layout(
title='Pareto curve',
xaxis_title='Cost',
yaxis_title='- Tastiness',
width=size, height=size,
)
fig.update_xaxes(range=[0, 11])
fig.show()
So, what does this line means? For every point in this line, there is no other feasible point that simultaneously reduces cost and -tastiness at the same time, meaning, choosing a point in this line means you had made a trade-off between cost and tastiness.
As one outstanding professor once said during my undergraduate studies "Life is a small blanket: you cover your head and you end up uncovering your feet, you cover your feet and you end up uncovering your head". Meaning, on multi-objective optimization problems in real life the objectives always end up competing, there is no "solution to rule them all", and you will often need to choose a trade-off between objectives.
Does the Pareto curve only exist in 2D? No, but for visualizing it we usually stop in 3D. As an example, let's add the minimization of calories as a third objective:
$$ \begin{equation} \tag{7} \begin{split} \text{minimize } \ & \mu_0 c^tx - \mu_1\tau^tx + \mu_2\sigma^tx\\ \text{subject to } \ & Ax \succeq b \\ & Ax \preceq 1.2b \\ & x \succeq 0 \end{split} \end{equation} $$In Eq. 7 the variables $\mu_0$, $\mu_1$, and $\mu_2$ will assume many values between $[0,1]$ so we can build the Pareto surface.
import itertools
A = df_processed.iloc[:, 3:].T.values
b = daily_intake.daily_recommended_intake.values.copy()
b_no_calories = b.copy()
b_no_calories[0] = 0
calories = df_processed.calories_k.values
c = df_processed.USD.values
n_points = 25
mu0 = [0] + np.logspace(-5, 0, num=n_points-1).tolist()
mu1 = [0] + np.logspace(-5, 0, num=n_points-1).tolist()
mu2 = [0] + np.logspace(-5, 0, num=n_points-1).tolist()
weight_combinations = itertools.product(mu0, mu1, mu2)
x = cp.Variable(shape=len(calories))
constraints = [x >= 0, A@x >= b_no_calories, A@x <= 1.2*b]
def optimize_diet(mu0, mu1, mu2):
objective_function = cp.Minimize(mu0*c@x - mu1*tastiness@x + mu2*calories@x)
problem = cp.Problem(objective=objective_function, constraints=constraints)
problem.solve()
cost_ = c@x.value
tastiness_ = -tastiness@x.value
calories_ = calories@x.value
return cost_, tastiness_, calories_
from joblib import Parallel, delayed
results = Parallel(n_jobs=-1)(
delayed(optimize_diet)(mu0, mu1, mu2)
for mu0, mu1, mu2 in tqdm(list(weight_combinations))
)
cost_list = [tuple_[0] for tuple_ in results]
tastiness_list = [tuple_[1] for tuple_ in results]
calories_list = [tuple_[2] for tuple_ in results]
100%|██████████| 15625/15625 [00:21<00:00, 739.86it/s]
import scipy.interpolate as interp
mesh_n_points = 200
x, y = np.meshgrid(
np.linspace(min(cost_list), max(cost_list), num=mesh_n_points),
np.linspace(min(tastiness_list), max(tastiness_list), num=mesh_n_points),
)
z = interp.griddata(
points=(cost_list, tastiness_list),
values=calories_list,
xi=(x, y),
method='linear'
)
fig = go.Figure(data=go.Surface(
x=x,
y=y,
z=z,
opacity=0.95
))
fig.update_layout(
title="Pareto 3D curve",
height=1000,
scene = dict(
xaxis_title = "Cost (x)",
yaxis_title = "- Tastiness (y)",
zaxis_title = "Calories (z)",
)
)
fig.show()
Although this diet optimization problem is well known and used for academic purposes, my objective with this blog post is to show the reader that defining the objective and translating it mathematically are not always easy tasks (heard of perverse instantion?), and, as most problems in real life, we will have competing objectives, meaning, choosing a trade-off between the multiple objectives is added to this already complex task.
That's it! I hope this was informative and made you think about properly setting objectives in your optimization problems. 💪